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Abstract 

Pulsed combustion is receiving renewed interest as a potential route to higher performance in air 
breathing propulsion systems. Pulsejets offer a simple experimental device with which to study unsteady 
combustion phenomena and validate simulations. Previous computational fluid dynamic (CFD) 
simulation work focused primarily on the pulsejet combustion and exhaust processes. This paper 
describes a new inlet sub-model which simulates the fluidic and mechanical operation of a valved pulsejet 
head. The governing equations for this sub-model are described. Sub-model validation is provided 
through comparisons of simulated and experimentally measured reed valve motion, and time averaged 
inlet mass flow rate. The updated pulsejet simulation, with the inlet sub-model implemented, is validated 
through comparison with experimentally measured combustion chamber pressure, inlet mass flow rate, 
operational frequency, and thrust. Additionally, the simulated pulsejet exhaust flowfield, which is 
dominated by a starting vortex ring, is compared with particle imaging velocimetry (PIV) measurements 
on the bases of velocity, vorticity, and vortex location. The results show good agreement between 
simulated and experimental data. The inlet sub-model is shown to be critical for the successful modeling 
of pulsejet operation. This sub-model correctly predicts both the inlet mass flow rate and its phase 
relationship with the combustion chamber pressure. As a result, the predicted pulsejet thrust agrees very 
well with experimental data. 


Nomenclature 


A e 

Exit plane area 

Ain 

Inlet cross sectional area 

A v 

Reed valve flow area 

A 

aa vjnin 

Minimal reed valve flow area 

gc 

Newton constant 

K 

Reed valve spring constant 

Lin 

Eength of the upstream end of the inlet 

m v 

Reed valve mass 

m 

Inlet mass flow rate 

m e 

Exit plane mass flow rate 

p 

1 cc 

Combustion chamber pressure 

Pin 

Pressure at the upstream end of the inlet 

Po 

Ambient total pressure 

Pup 

Pressure just upstream of the reed valve 
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V e Exit plane average velocity 

V in Inlet velocity 

V v Reed valve velocity 

x v Reed valve position 

a Proportionality constant relating open reed valve area to its position 

p density 

x Thrust 


Introduction 

The pulsejet is of interest both due to the complexity of the fundamentally unsteady, reactive 
flowfield and due to its potential applications to create low-cost, scalable, high thrust-to-weight ratio 
propulsion devices (refs. 1 and 2). There is also renewed interest in the application of pulsed combustion 
devices for achieving pressure-gain heat addition in gas turbine systems (ref. 3). 

Detailed descriptions of the mechanically-valved pulsejet operational cycle can be found in the 
literature (refs. 4 and 5); however, it may briefly be described as follows. Referring to figure 1, the cycle 
may be divided into three phases 

(1) Combustion — Reaction of an air and fuel mixture within the combustion chamber commences. 

The pressure begins to rise as a result of confinement of the flow. The pressure rise causes the inlet valves 
to close, preventing backflow. The reaction accelerates as the pressure and temperature rise; this, in turn, 
accelerates the pressure and temperature rise. 

(2) Expansion — The hot, high-pressure gases in the combustion chamber expand, forcing flow 
through the exhaust. 

(3) Ingestion — The momentum of the exhaust gases causes the combustion-chamber pressure to 
drop below the ambient value. This allows the inlet valves to open and a fresh charge of air to enter 
(mixed with fuel). Eventually, the exiting exhaust flow reverses and mixes with the fresh charge. This 
initiates a new reaction, and the cycle begins again. 

The design and optimization of any pulsed combustion system will rely on accurate flowfield 
simulations. However, the flowfield, both internal and external, is very complex. For example, the 
exhaust flow contains a strong starting vortex common to all unsteady, thrust producing flows. The 
starting vortex has been extensively studied both experimentally and computationally (refs. 6 to 13), in a 
variety of unsteady devices. However, in most of that work, the flowfields are not nearly as complex as 
that of the pulsejet. That is to say, they do not contain comparable turbulence levels, driver and driven 
enthalpy ratios, vorticity levels, or accelerations. Other challenges include predicting the periodic 
combustion process, the associated pressure rise in the combustion chamber, and the complex, coupled 
gasdynamic processes in the tailpiece that result in resonant operation. The inlet mass flow, which is a 
function of the spring force of the reed valve, the combustion chamber pressure and the upstream flow 
properties, needs to be modeled precisely to provide appropriate response to the internal flowfield 
variations between cycles. In other words, simulations should be able to model the self-aspirating 
phenomena of pulsejets. There have been limited numerical studies of pulsejet flows (e.g., refs. 14 to 17), 
and some experimental flowfield measurements (beyond just thrust) (ref. 18). However, the two have 
rarely been directly compared, and without such comparisons, or validation efforts, conclusions drawn 
from the computational work are less certain. This paper represents a first step toward that end. 
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Using an experimentally guided, computational fluid dynamic approach, previous work (refs. 17 to 
19) established simulation models and studied the operation of various pulsejet configurations. Numerical 
simulations were performed utilizing CFX, a commercially available, 3-D, compressible, viscous, 
computational fluid dynamics (CFD) code. In order to capture wave motions and their interactions and to 
minimize numerical diffusion in the simulations, second order advection schemes were always used to 
calculate the advection terms in the discrete finite volume equations. The CFD mesh and time step were 
optimized to ensure mesh-independent results. These simulations provided physical insights into the 
pulsejets operation. Pressure, mean temperature, exit velocity, static thrust, and operational frequency 
were compared between simulation results and limited experimental data. Although the results showed 
good agreement in pressure, temperature, exhaust plane velocity, and frequency, simulated total thrust 
and inlet mass flow rate were less than those obtained in experiments. This was due to the fact that the 
pulsejet valve head, which essentially consisted of a short duct and a reed valve, was not fully modeled. 
Instead the flow through the valve was modeled using an inlet boundary condition (ref. 17) that defined 
the inlet velocity (at standard temperature and pressure) as a linear function of the pressure difference 
between the ambient pressure and the combustion chamber pressure. This function, based on 
experimental measurements, is given below: 


V in = 


3(P 0 -P cc )/mO;P cc <P 0 

0;P CC >p 0 


( 1 ) 


where P cc is the pressure immediately after the reed valves and P 0 is the ambient pressure. Clearly this 
simple function could not model the inertia effect of the upstream flow in the valve head and valve 
motions, resulting in a poor prediction of the inlet mass flow rate and its phase relation to the time 
varying combustion chamber pressure. 

The objective of this paper is to improve the simulation model through a modified sub-model of the 
inlet and reed valves. Similar to the previous work, pressure and velocity data will be compared to 
experimental data. However, this paper will concentrate on the flowfield behind or downstream of the 
pulsejet exit plane. This is the region where the starting vortex is formed, and arguably represents the 
most intriguing and challenging region of the flow. Furthermore, the emitted vortex has been strongly tied 
to the very high levels of entrainment, mixing, and thrust augmentation observed when pulsejets are used 
in conjunction with ejectors (refs. 20 and 21). Maximizing this augmentation represents enabling 
technology for many potential pulsejet-based applications. It is therefore critical to validate the developed 
simulation code in this flow region. 

Furthermore, a detailed set of experimentally measured, velocity related quantities have been 
obtained in this region and can serve as a set of valuable validation data; particularly since the simulated 
pulsejet was identical to the one used for the velocity measurements. The experimental data, obtained 
using particle imaging velocimetry (PIV), were collected in the same manner described in reference 13 
(and will therefore not be presented here). In fact, the same equipment and laboratory setting were used. 
However, for the reference 13 measurements, an ejector was present in the exhaust region. For the present 
work, the ejector was removed. 


Simulation Model 

Simulation Model and Grid 

Unless specifically stated in this paper the simulation model shares similar features as those in 
previous work (refs. 17 to 19). Figure 2 shows the dimensions of the computational domain used in 
simulating the pulsejet and external flowfield. These dimensions match those of the experimental device 
to which comparison will be made. The model assumes that the flow is 2-D axis-symmetric and only 1/90 
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Figure 2. — Dimensions of the simulated pulsejet flowfield (mm). 



Figure 3. — Hybrid mesh for the pulsejet simulation. 


in the 0-direction of the whole domain is modeled. The turbulent flow and premixed propane combustion 
are modeled with the Shear-Stress-Transport (SST) model and Eddy Dissipation Model, respectively 1 . 
The k - co based SST model is used to give highly accurate predictions of the onset and the amount of flow 
separation under adverse pressure gradients by the inclusion of transport effects into the eddy-viscosity 
formulation. The result is a major improvement in terms of predicting flow separation, which is expected 
to occur in the flared region of the pulsejet tailpiece. The effectiveness of this approach has been 
demonstrated in a large number of validation studies (ref. 22). The viscous effect in the boundary layer 
region is modeled by the CFX automatic near-wall treatment which automatically switches from wall- 
functions to a low-Reynolds near wall formulation as the mesh is refined. However, to take advantage of 
the reduction in errors offered by this low-Reynolds formulation, at least 10 nodes must be put in the 
boundary layer, causing a considerable increase of mesh size. A hybrid high quality mesh is created to 
ensure fine hex-mesh in the near wall region and exit plane while relatively coarse tet-mesh in part of the 
external flowfield, as shown in figure 3. 


! The experimental pulsejet to which comparison will be made was actually liquid- fueled, with gasoline. However, 
accurately modeling the complexities of liquid droplet sprays, vaporization, etc. was considered beyond the scope of 
this work. 
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There are total 13766 elements, of which 362 are tet-meshes. The computations are performed on the 
North Carolina State IBM Blade center utilizing a single 3.0 GHz Intel Xeon processor. With a 1.0 psec. 
time step, the computation time for one cycle is about 24 hr. 


Inlet Sub-model 

The computational simulation includes a modified sub-model of the inlet and reed valve section that 
is more realistic than previously used (ref. 17). This sub-model is implemented as a boundary condition 
controlling the mass flow rate of the premixed reactants into the combustion chamber. As shown in 
figure 4, the mass flow rate is a function of pressure in the valve head and combustion chamber pressure. 
The valve motion is a function of spring force and pressure forces acting on both sides of the valve. These 
two processes are coupled in the sub-model. The inlet flow and reed valve dynamics are modeled by 
solving three ordinary differential equations that describe their motions. 


- 12^-g ( p . -p ) 
dt L,„ 8A up) 


( 2 ) 


dV v 

dt 



(3) 


The quantities P up , 


dx v 

dt 


■ = K 


and A v are found from the following relations. 


P,+- 


m 


P = 

up 


12g e 2p 

•2 


1 1 


P -- 


A" An 


, m < 0.0 


m > 0.0 


12g,2p,(,.- 

A v = max (A v min ,a x v ) 


(4) 

(5) 

( 6 ) 


Equation (2) describes the mass flow rate and equations (3) and (4) describe the valve motion. The 
combustion chamber pressure, P cc , is the force function calculated by the code. At each timestep, the sub- 
model calculates the new mass flow rate based on the current combustion chamber pressure and returns 
this value to the code, where inlet mass flow rate is specified. Equations (2) through (4) are integrated 
using a simple, first-order Euler method, with a time step equal to that of the pulsejet simulation. 
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A/WM 

k i 


K v 


m 


Figure 4. — Valve head model schematic. 
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Figure 5 shows the sub-model response to a sinusoidal variation of the downstream combustion 
chamber pressure (representing resonant pulsejet operation). The simulation results show that the 
maximum mass flow rate across the valve does not occur at the lowest combustion chamber pressure, as 
might be expected and implied in the previous work. Furthermore, flow into the inlet does not stop when 
the combustion chamber pressure exceeds the ambient pressure (101325 Pa). This sub-model has been 
partially validated by comparing simulated and experimentally measured reed valve positions and 
combustion chamber pressures during limit cycle operation, as explained next. 

Figure 6 shows the experimental data of combustion chamber pressure and valve position over the 
course of one cycle. This trace should be compared to figure 5(b). The pressure transducer installation is 
described in reference 20. The valve position is measured using specially fabricated fiber-optic sensor 
which is described in the Appendix. The sensor is not calibrated, since the design is crude, and accuracy 
is not expected. The intention of the measurement is to obtain information on when the valve opens and 
closes relative to the combustion chamber pressure. The precise position is not particularly important. As 
such, the position data is shown simply as sensor voltage subtracted from the voltage reading when the 
valve is in the closed position (sensor voltage drops when the valve opens). The experimental data of 
figure 6 looks quite similar to that of figure 5(b). In particular, the valve is still open after the combustion 
chamber pressure exceeds the ambient pressure, and it remains closed even after the combustion chamber 
pressure drops below ambient. Thus, the valve head sub-model and experimental similarities, while not 
conclusive, do suggest that the sub-model approach used here is valid. 




, A 

(a) (b) 


Figure 5. — Inlet and valve model limit cycle responses to a sinusoidal approximation of the pulsejet 
combustion chamber pressure: (a) upstream pressures and mass flow, (b) reed valve position. 
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Figure 6. — Measured pulsejet combustion chamber pressure and 
reed valve position. 
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Figure 7. — Simulated pulsejet inlet mass flow rate and 
combustion chamber pressure. 

Figure 7 presents the pulsejet simulation (with the valve head sub-model implemented) combustion 
chamber pressure and inlet mass flow rate over the course of one limit cycle. In the simulation model, 
negative flow (flow going from the combustion chamber to the valve head) is not allowed due to the 
numerical instability caused by the negative flow. Comparing figure 7 with figure 5(a), it is evident that 
the inlet sub-model is correctly implemented. The peak mass flow rate does not occur at the lowest 
combustion chamber pressure, but rather about 0.0005 s later. More importantly, the cycle-averaged inlet 
mass flow rate increases from 0.012 kg/s in previous work to 0.021 kg/s. The experimentally measured 
mass flow rate is 0.026kg/s (ref. 20), providing further indication that the valve head sub-model is correct. 
This improvement, increasing the inlet mass flux by 75 percent, is due solely to the utilization of the inlet 
sub-model and is critical to accuracy of the final simulation results because with increased mass flow rate 
of reactants, the total thrust increases. This is a desirable effect because previously calculated thrust was 
smaller than the experimentally measured thrust. Additionally, the correct phase relationship between 
inlet mass flow rate and combustion chamber pressure is expected to improve prediction of flow 
dynamics in and out of the pulsejet, which will be evident when examining the exhaust flow later. 


Simulated and Experimental Flow Comparisons and Discussion 

As stated in the introduction, the objective of this paper is to validate simulation results by comparing 
them with experimental data from the pulsejet. Such comparisons are made in this section and include 
global quantities, such as cycle-averaged thrust and mass flow, as well as time-dependant, local quantities 
such as combustion chamber pressure, and exhaust field velocities. 

The experimentally measured thrust is 4.3 lbf, or 19.13 N, with a standard deviation of 1.4 percent 
(ref. 20). The simulated result is shown in figure 8. Thrust is calculated as the integral of the momentum 
force at the pulsejet exit plane, i.e. 


t = J m e V e dA e (8) 

For both the previous and current results, the thrust profiles for three cycles are shown in this figure; 
cycle-averaged thrusts are computed to be 19.46 and 16.44 N for current and previous results, 
respectively. Compared to measurements, the error is 1.7 percent for current result and 14 percent for 
previous result. 
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Figure 8. — Simulation results for the pulsejet thrust. 
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Figure 9 shows the simulated and experimentally measured centerline exhaust velocities, at 7 and 
210 mm downstream of the pulsejet exit plane, over the course of one cycle. The experimental data give 
the average value over an approximately 6 mm span across the centerline, while the simulation data give 
the exit velocity exactly at the centerline. Combustion chamber pressure is also included in the plot in 
order to examine the pressure- velocity phase relationship. Time zero represents the time when 
combustion chamber pressure reaches its peak value. In experiments, the PIV data are taken at 125 psec. 
intervals commencing approximately 300 psec after time zero. This 300 psec delay is therefore applied to 
all the simulation results in this paper. It is evident that simulation agrees well with measurements in 
terms of pressure, velocity, and their phase difference. The simulated velocity at 210 mm downstream of 
the exit is higher than the measured one. This may be caused by the fact that flow velocity diffuses as it 
travels downstream. Thus, the average value over several points by the centerline may be smaller than the 
centerline velocity. 

Figure 10 compares computed and measured vorticity contours at three different instants during the 
pulsejet cycle. In each plot, the upper half is the measured PIV data and the lower half is the simulated 
data. The “+” sign in the PIV data represents the location of the highest vorticity value observed in the 
vortex region. This value, shown on the right hand side of the PIV contour, is used to scale the color 
spectrum shown, thereby accentuating the vortex. It can be seen from these plots that the simulated vortex 
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has the same convective speed as the one in measurements. The vorticity levels also agree well between 
simulation result and experimental data. Figure 10(a) shows that the starting vortex ring has just formed 
from the pulsejet exit shear layer and is still attached to the trailing jet. Possibly owing to this attachment, 
the vortex in figure 10(a) has the highest vorticity value in the vortex core. As the vortex travels 
downstream, it starts to detach from the trailing jet. As a result of this, and dissipative forces, the vorticity 
level of the vortex ring begins to decrease, as shown in figure 10(b). Later, in figure 10(c), the vortex is 
completely detached from the trailing jet, and the vorticity level continues to decrease. 
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Figure 10. — Experimentally measured and simulated vorticity 
contours at (a) T = 0.000375 s, (b) T = 0.00075 s, and 
(c)T = 0.0015 s. 
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Figure 11. — Measured and simulated vortex location and 
combustion chamber pressure versus time. 


To further verify the simulation results, computational and experimental data on vortex location and 
combustion chamber pressure are compared in figure 1 1 . The location of the vortex is represented by the 
axial distance from the exit plane to the point in the vortex with the highest vorticity value. From time 
zero, the vortex locations at each 125 psec. interval are compared between simulation and experiment. 

The comparison shows a good agreement with small but consistent differences after 0.00125 s. One 
possible reason for this difference is that the vortex has traveled to the edge of the simulation domain, 
where a constant pressure of 1.01 kPa was specified. This plot also reveals that for the first 0.0015 sec., 
although the velocity at the centerline changes from about 500 to 200 m/s (fig. 9), the vortex ring travels 
at a relatively constant speed of about 135 m/s according to the data. 

Conclusion 

This paper provides direct comparisons between experimental measurements and simulation results of 
mechanically valved pulsejet. The simulation contains a new inlet sub-model which simulates the valve 
head dynamics and is integrated into the code. As a result, this new model is able to correctly predict the 
total thrust and time-averaged mass flow rate, an improvement over the previous model. 

Besides thrust, the simulation model is fully validated by comparisons with experimental data on 
combustion chamber pressure, exit flow velocity, pressure-velocity phase relationship, vortex ring 
location and vorticity level. The results obtained indicate that the peak inlet mass flow rate occurs about 
l/9th of a cycle after the minimum in combustion chamber pressure. The inertia effect of the flow inside 
the valve head plays a key role in the inlet mass flow rate modeling. Without inertia, the model tends to 
underpredict the inlet mass flow rate, by as much as 75 percent. 

The comparison of the simulated pulsejet exit flowfield with that measured experimentally also 
verifies simulation results. It may well be possible therefore, to utilize the developed simulation to gain 
insights into the role of the emitted vortex in achieving substantial thrust augmentation and mixing, when 
an ejector is present behind the pulsejet. 
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Appendix — Optical Sensor Detail 

The optical position sensor constructed for this work is shown installed in the composite photograph 
of figure 12. An infrared LED transmits light down one of the fiber-optic rods shown. This light is 
emitted from the end of the rod, which is placed in close proximity to the face of the reed valve, in the 
closed position. The emitted light is therefore reflected from the valve face, and some of it is transmitted 
through the second fiber optic rod, which terminates into an infrared phototransistor. The light is 
sufficient to allow the transistor to pass current which, passing also through a resistor, registers a high 
voltage signal. When the reed valve opens, less of the reflected LED light reaches the transistor, causing 
the signal voltage to drop. 



Figure 12. — Optical reed valve position sensor. 
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